##############################################################################
### replicating the R-based figures from:
### Bräuninger, T., Däubler, T., Huber, R., & Rudolph, L. 
### How Open Lists Undermine the Electoral Support of Cohesive Parties
### forthcoming in British Journal of Political Science
##############################################################################

library(readstata13)
library(dplyr)
library(ggplot2)

########################
### Figures 2 and 3  ###
########################

library(basicspace)
library(grid)
library(gridExtra)

# German Longitudinal Election Study (2019). Post-election Cross Section (GLES 2017). GESIS Data Archive, Cologne. 
# ZA6801 Data file Version 4.0.1
# available at : https://doi.org/10.4232/1.13235

# d <- read.dta13("../GLES_analysis/ZA6801_de_v4-0-1.dta", convert.factors = FALSE)
                
d <- select(d, c(bula, q76a:q76g, q79)) #  state/Bundesland; party perceptions; self-placement
d[,2:9] <- apply(d[,2:9],2, function(x) ifelse(x < 0, NA, x)) # recode negative values = missings
names(d)[2:9] <- c(paste0("ptypercep",1:7), "self")

##  perceptions of parties and self placement

# for CDU/CSU: use CSU if Bavaria, CDU otherwise
dpercselfBav <- apply(select(d,c(ptypercep1:ptypercep7, self)),2,as.numeric)

dpercselfBav[d$bula==9,1] <- dpercselfBav[d$bula==9,2]
dpercselfBav <- dpercselfBav[,c(1,3:8)]
iscompleteBav <- rowSums(!is.na(dpercselfBav))==7
dim(dpercselfBav)

colnames(dpercselfBav) <- c("CDU/CSU", "SPD" , "Left" , "Greens" ,"FDP" , "AfD","self")

# Share of missing in party placement 
dim(dpercselfBav)
summary(dpercselfBav)
mis <- apply(dpercselfBav[,1:6],2, function(x) sum(is.na(x)/length(x)) )
names(mis) <-  c("CDU/CSU","SPD" , "Left" , "Greens" ,"FDP" , "AfD")
mis


### Aldrich-McKelvey scaling

amBav <- aldmck(dpercselfBav, respondent = 7,  polarity= 4, verbose=T)
sdPlac <- apply(dpercselfBav[,1:6], 1, sd, na.rm=TRUE) # no respondent params for people placing all parties at same point

dsc <- as.data.frame(dpercselfBav[,1:6] * amBav$respondents$weight + amBav$respondents$intercept)
dsc <- dsc[amBav$respondents$weight > 0,]

dsc$id <- 1:nrow(dsc)
names(dsc) <- c(paste0("var_",1:6),"id")
dsc <- reshape(dsc, idvar="id", varying=1:6, sep="_", direction="long")
dsc$Party <- as.factor(dsc$time)
levels(dsc$Party) <- c("CDU/CSU","SPD" ,"Left" , "Greens" ,"FDP" , "AfD")
dsc$Party <-  factor(dsc$Party,levels(dsc$Party)[c(3,4,2,5,1,6)])

# this is for respondents placing all parties and themselves (and sd of placements != 0)
table(iscompleteBav & sdPlac > 0 & amBav$respondents$weight > 0) 
nrow(dsc[!is.na(dsc$var),])/6

pmain <-ggplot(dsc[!is.na(dsc$var),],aes(Party,var)) +
  geom_boxplot() + coord_flip() +
  labs(y="Scaled party placements") +
  annotate("text", x=1:6, y = -1, label=paste0("(",100*round(mis[c(3,4,2,5,1,6)],3),")"), size=4.5) +
  theme_bw() + theme(axis.text=element_text(size=12), # text of labels
                     axis.title=element_text(size=14)) 
pmain 
ggsave(file="figures/Fig2.png", width=9, height=6)

den <- density(amBav$respondents$idealpt[amBav$respondents$weight >0], n =2^13, na.rm=T)
#plot(den)
dp <- data.frame(x=den$x, d=den$y, cd=cumsum(den$y)/max(cumsum(den$y)))
tmp <- plyr::laply(1:6, function(n) which(abs(dp$x-amBav$stimuli[n])==min(abs(dp$x-amBav$stimuli[n]))) )
dp2 <- dp[tmp,]
dp2$party <- colnames(dpercselfBav)[1:6]

p1 <- ggplot(dp, aes(x,d)) +
  geom_line() +
  geom_segment(data=dp2,aes(xend=x,yend=0),linetype="dotted") +
  geom_text(data=dp2,aes(x=x+c(.05,0,.04,-.06,-.05,0),y=-.01,label=party), angle=90,hjust=1, size=4) +
  xlim(-2,2) + ylim(-.11,.805) + 
  labs(x="Scaled immigration positions", y="",title="Density") + theme_bw() + 
  theme(axis.text=element_text(size=12), # text of labels
        axis.title=element_text(size=14)) 

p2 <- ggplot(dp, aes(x,cd)) +
  geom_line() +
  geom_segment(data=dp2,aes(xend=x,yend=0),linetype="dotted") +
  geom_segment(data=dp2,aes(xend=-1.8,yend=cd),linetype="dotted") +
  geom_text(data=dp2,aes(x=-1.8,y=cd+c(.03,.03,.03,-.03,-.03,.03),label=party), hjust=0, size=4) +
  xlim(-1.8,1.8) + labs(x="Scaled immigration positions",y="", title="Cumulative density") + theme_bw() +
  theme(axis.text=element_text(size=12), # text of labels
        axis.title=element_text(size=14)) 

png(file="figures/Fig3.png", width=12, height=6, units="in", res=300)
grid.arrange(arrangeGrob(p1, p2, ncol = 2), 
             nrow = 1) 
dev.off()

rm(d)


##############################
### Figure 4 and Figure A2 ###
##############################

d <- read.dta13("data.dta")

d$revimmig <- -d$immigration # reverse scale

sm <- select(d,starts_with("SM"))
stopifnot(!is.na(sm)) # no missings here

rowmax <- apply(sm,1,max)

smeqmax <- apply(sm,2,function(x) as.numeric(x==rowmax)) # does entry equal max? 0/1
natmax <- rowSums(smeqmax) # how many at max?
smnomax <- sm
smnomax[smeqmax==1] <- NA
rowsecond <- apply(smnomax,1,max,na.rm=TRUE)
smtop2 <- apply(sm,2,function(x) as.numeric(x==rowmax | (x==rowsecond & x >= rowmax-1  )))
smtop2[natmax > 1,] <- smeqmax[natmax > 1,] # do not consider rowsecond if more than one party at max 
# smtop2 = is in top2 0/1 (can be true for more than 2 parties if several ranked at max)

d$indiffCDUAFD <- (smeqmax[,3]==1 & sm[,1] >= rowmax-1) | (smeqmax[,1]==1 & sm[,3] >= rowmax-1)
prop.table(table(d$indiffCDUAFD))
table(d$indiffCDUAFD[d$treatmentR1 != "OpenInfo+"])

d$indiffSPDGRU <- (smeqmax[,4]==1 & sm[,2] >= rowmax-1) | (smeqmax[,2]==1 & sm[,4] >= rowmax-1)
prop.table(table(d$indiffSPDGRU))
table(d$indiffSPDGRU[d$treatmentR1 != "OpenInfo+"])

dg <- summarize(group_by(d, revimmig),
                shIndiffCDUAFD = sum(indiffCDUAFD),
                shIndiffSPDGRU = sum(indiffSPDGRU),
                N = n()) # note that for now shIndiff is the absolute number
# calculate CIs using prop.test
dg$lowerCDUAFD <- NA
dg$upperCDUAFD <- NA
dg$lowerSPDGRU <- NA
dg$upperSPDGRU <- NA


for(i in -10:0){
  tmp <- prop.test(dg$shIndiffCDUAFD[dg$revimmig==i], dg$N[dg$revimmig==i])
  dg$lowerCDUAFD[dg$revimmig==i] <- tmp$conf.int[1]
  dg$upperCDUAFD[dg$revimmig==i] <- tmp$conf.int[2]
  tmp <- prop.test(dg$shIndiffSPDGRU[dg$revimmig==i], dg$N[dg$revimmig==i])
  dg$lowerSPDGRU[dg$revimmig==i] <- tmp$conf.int[1]
  dg$upperSPDGRU[dg$revimmig==i] <- tmp$conf.int[2]
}
dg$rowshIndiffCDUAFD <- dg$shIndiffCDUAFD/dg$N
dg$rowshIndiffSPDGRU <- dg$shIndiffSPDGRU/dg$N


ggplot(dg, aes(revimmig,rowshIndiffCDUAFD)) +
  geom_point(size=3) + geom_line(size=1) +
  geom_errorbar(aes(ymin=lowerCDUAFD, ymax=upperCDUAFD), width=.2) +
  scale_x_continuous(breaks=seq(-10,0,2), labels=(seq(0,10,2))) +
  scale_y_continuous(limits=c(0,.20), breaks=seq(0,.18,.06)) +
  labs(x="Position on immigration \n (0=facilitate immigration, 10=restrict immigration)",
       y="Share among all respondents at position" ) +  
  theme_bw() + theme(axis.text=element_text(size=12), # text of labels
                 axis.title=element_text(size=14)) 
ggsave(file="figures/Fig4.png", width=9, height=6)

ggplot(dg, aes(revimmig,rowshIndiffSPDGRU)) +
  geom_point(size=3) + geom_line(size=1) +
  geom_errorbar(aes(ymin=lowerSPDGRU, ymax=upperSPDGRU), width=.2) +
  scale_x_continuous(breaks=seq(-10,0,2), labels=(seq(0,10,2))) +
  scale_y_continuous(limits=c(0,.36), breaks=seq(0,.35,.05)) +
  labs(x="Position on immigration \n (0=facilitate immigration, 10=restrict immigration)",
       y="Share among all respondents at position" ) + 
  theme_bw() + theme(axis.text=element_text(size=12), # text of labels
                     axis.title=element_text(size=14)) 
ggsave(file="figures/FigA2.png", width=9, height=6)

########################
### Figures 7 and A6 ###
########################

library(gtable)

d <- read.dta13("data.dta")

### select subset - here Info (narrow) groups only
d <- d[d$treatmentR1 %in% c("ClosedInfo", "OpenInfo"),  ]

d$open <- ifelse(d$treatmentR1 ==  "OpenInfo", 1, 0)

d$afdR1 <- ifelse(d$votechoiceR1=="AfD", 1, 0)
d$gruR1 <- ifelse(d$votechoiceR1=="GRUENE", 1, 0)

# reverse immig. scale
tmp <- data.frame(immigration=0:10, tmp=10:0)
d <- merge(tmp,d)
d$immigration <- d$tmp

sm <- select(d,starts_with("SM"))
rowmax_exclAfD <- apply(sm[-c(3)],1,max)
d$indiffAfD <- ifelse(abs(sm$SM_AfD - rowmax_exclAfD)<2, "Indifferent", "Preference")
rowmax_exclGru <- apply(sm[-c(4)],1,max)
d$indiffGru <- ifelse(abs(sm$SM_Gru - rowmax_exclGru)<2, "Indifferent", "Preference")

table(d$treatmentR1, d$immigration)

d$ruler <- ifelse(d$open==0, d$immigration-.25, d$immigration+.25)
d$ballotXafd <- 1
d$ballotXafd[d$open == 0 & d$afdR1 == 1] <- 2
d$ballotXafd[d$open == 1 & d$afdR1 == 0] <- 3
d$ballotXafd[d$open == 1 & d$afdR1 == 1] <- 4

## data: A[fd]I[ndifferent],G[reens]I[ndifferent]
dAI <- d[d$indiffAfD == "Indifferent",]
dGI <- d[d$indiffGru == "Indifferent",]


mytest <- function(imm,dat, greens=FALSE, OminusC=FALSE){
  
  if(greens==TRUE) choice <- dat$gruR1
  else choice <- dat$afdR1
  tab <- table(dat$open[dat$immigration==imm], 1-choice[dat$immigration==imm])
  
  if(nrow(tab) != 2 | ncol(tab) != 2)  {
    out <- rep(NA,9)
  }  else {
    tmp <- prop.test(tab)
    out <- c(tmp$estimate[1], tmp$estimate[2], tmp$estimate[1]-tmp$estimate[2],tmp$conf.int[1], tmp$conf.int[2])
    absnum <- c(sum(choice[dat$open == 0 & dat$immigration==imm]), sum(1-choice[dat$open == 0 &  dat$immigration==imm]),
                sum(choice[dat$open == 1 & dat$immigration==imm]), sum(1-choice[dat$open == 1 &  dat$immigration==imm]))
    if(OminusC==T) out <- c(tmp$estimate[1], tmp$estimate[2], tmp$estimate[2]-tmp$estimate[1],-tmp$conf.int[1], -tmp$conf.int[2])
    out <- c(absnum, out)
  }
  names(out) <- c("absChom","absCnhom","absOhom","absOnhom","closed","open","diff","lower","upper")
  return(out)
}

testresAI <- plyr::ldply(0:10,  mytest, greens=FALSE, dat=dAI, OminusC=T)
testresGI <- plyr::ldply(0:10,  mytest, greens=TRUE, dat=dGI, OminusC=T)

testresAI$immigration <- 0:10
testresAI$N <- rowSums(testresAI[,1:4])
testresAI[2,c(1:4,6,11)] <- c(0,0,1,1,.5,2) # fix entries due to incomplete table
testresAI$open[testresAI$open==0] <- .005 # add a bit so that this is plotted

testresGI$immigration <- 0:10
testresGI$N <- rowSums(testresGI[,1:4])

mybreaks <- c(0:10 - .25, 0:10 + .25, 0:10)
mylabels <- c(rep("C",11),rep("O",11),as.character(0:10))

myplots <- function(testres, plab, plotname){
  
  tmp2 <- testres[,c("closed","open","immigration","N" )]
  names(tmp2) <- c("val_0", "val_1","immigration","N" )
  tmp2 <- reshape(tmp2, varying=1:2, sep="_", idvar="immigration", timevar="open",direction="long")
  tmp2$ruler <- ifelse(tmp2$open==0, tmp2$immigration-.25, tmp2$immigration+.25)
  tmp2$List <- factor(1-tmp2$open)
  levels(tmp2$List) <- c("open","closed")
  
  pbarrel <- 
    ggplot(tmp2) +
    geom_col(aes(x = immigration, y = val, fill = List), position = "dodge", width=.7) +
    geom_line(aes(x=immigration, y=N/400),color="black",size=1.2) +
    scale_x_continuous(breaks=0:10, limits=c(-.5,10.5)) + 
    scale_y_continuous(limits=c(0,.8), sec.axis = sec_axis(~.*400, name = "N across list types (line)")) + 
    scale_fill_manual(values=c("grey35","grey")) +
    labs(x = "",y=paste0("Vote proportion for ",plab," (bars)")) +
    theme_bw()   +  theme(legend.position = "bottom", legend.box.margin=margin(-25,0,0,0)
    )
  
  pmain <- 
    ggplot(testres, aes(x=immigration, y=diff, ymin=lower, ymax=upper)) +
    geom_point() + geom_errorbar(width=.2) + 
    geom_line() + # colour="blue"
    geom_hline(yintercept=0, linetype="dashed") +
    scale_x_continuous(breaks=0:10, limits=c(-.5,10.5)) +
    labs(x="Position on Immigration\n (0=facilitate immigration, 10=restrict immigration)",
         y="Difference in proportions:\n open - closed") +
    theme_bw() +
    theme(panel.grid.minor=element_blank()) 
  
  g1 <- ggplotGrob(pbarrel)
  g2 <- ggplotGrob(pmain)
  
  g <- rbind(g1, g2, size = "first")
  g$widths <- unit.pmax(g1$widths, g2$widths)
  grid.newpage()
  png(file=paste0("figures/",plotname,".png"), width=12, height=16, units = "cm", res=600)
  grid.draw(g)
  dev.off()
}

myplots(testresAI, plotname="Fig7",plab="AfD")
myplots(testresGI, plotname="FigA6",plab="Greens")

##################
### Figure A7  ###
##################

# Use second wave of Comparative Candidate Survey to look at immigration position of candidates

# data is available here:
# https://forsbase.unil.ch/project/study-public-overview/17220/0/

#d <- read.dta13("886_CCS_Data_Wave2_v4.0.0.dta")

# attention: GER and ICE have more than one year/survey in here
table(d$T1)
library(stringr)
d$country <- str_replace_all(as.character(d$T1)," \\d{4}","")
d$year <- str_extract_all(as.character(d$T1),"\\d{4}","")

d <- d[!(d$country %in% c("Albania","Australia","Canada","Chile","Montenegro")),]
table(d$country, d$T8)
d <- d[!(!is.na(d$T8) & d$T8 == "constituency candidate"),] # drop SMD-only candidates in mixed-member system

# "Immigrants should be required to adapt to the customs of [country]"
table(d$C2a)
d$adapt <- as.numeric(d$C2a)-3
table(d$adapt)

# Is candidate voting possible (flexible or open)?
d$candvot <- 1 # includes non-list: Ireland
d$candvot[d$country %in% c("Norway","Portugal","Romania","Spain", "Germany","Hungary")] <- 0
# includes closed-list mixed (GER, HUN)

table(d$candvot, d$year)

# aggregate
library(DescTools)
da <- summarize(group_by(d, candvot, country, year, A1),
                extAdapt=mean(adapt==1,na.rm=T),  # share in the most rightist category
                entAdapt=Entropy(table(adapt)),
                NAdapt=sum(!is.na(adapt)),
)
da <- da[da$NAdapt > 9,] # require at least 10 answers per party
cor(da[,5:6])

da$CandidateVoting <- as.factor(da$candvot)

g1 <- ggplot(da, aes(extAdapt, entAdapt, fill=CandidateVoting)) +
  geom_point(shape=21, colour="black") +
  scale_fill_manual(values=c("white","black")) +
  labs(x="Proportion in most rightist of 5 categories", y="Entropy") +
  theme_bw() + theme(aspect.ratio=1, legend.position="bottom")
g1
ggsave(file="figures/FigA7.png",width=24, height=10, units="cm", dpi=600)


